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Abstract 



'!;;j- ' The dynamics of the normal/superconducting interface in type-I supercon- 

ductors has recently been derived from the time-dependent Ginzburg-Landau 
theory of superconductivity. In a suitable limit these equations are mapped 
^ ' onto a "free-boundary" problem, in which the interfacial dynamics are deter- 

o' 

J-^ ' mined by the diffusion of magnetic flux in the normal phase. The magnetic 

L^ ■ field at the interface satisfies a modified Gibbs-Thomson boundary condition 

C^ ' which involves both the surface tension of the interface and a kinetic coeffi- 

cient for motion of the interface. In this paper we calculate the surface tension 
and kinetic coefficient numerically by solving the one dimensional equilibrium 
Ginzburg-Landau equations for a wide range of k values. We compare our 
numerical results to asymptotic expansions valid for k ^ 1, k ~ l/v2, and 

K ^ 1, in order to determine the accuracy of these expansions. 
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I. INTRODUCTION 

When a type-I superconductor in a magnetic field is subjected to either a sudden tem- 
perature or magnetic field quench which takes it from the normal phase into the Meissner 
phase, the approach to equilibrium will be determined by the rate at which superconducting 
islands are nucleated in the background normal phase, and the subsequent dynamics of the 
superconducting/normal interfaces. Refs. [|I| and suggested that the essential features of 
the interface motion could be understood in terms of a free boundary model for the magnetic 
field in the normal phase; this model is almost identical to a free boundary model which is 
used to study the growth of a solid into its supercooled liquid phase. The interface motion 
in the latter case is known to be unstable, and leads to highly ramified solidification pat- 
terns (dendrites, for instance; see Ref. for an overview). By analogy, the growth of the 
superconducting phase into the normal phase should be dynamically unstable. Numerical 
solutions of the time- dependent Ginzburg-Landau (TDGL) equations of superconductivity 
confirmed these expectations |l|J^. However, the precise connection between the TDGL 
equations and the free boundary model was not made. 

More recently, the free boundary model has been derivedirora the TDGL equations using 
the method of matched asymptotic expansions [^,Q . The free boundary model consists of a 
diffusion equation for the magnetic field h in the normal phase, 

dth = DHV\ (1.1) 

where Dh = l/A-na^'^^ is the diffusion constant for the magnetic flux, with a^^'^ the 
normal state conductivity; a continuity equation for the magnetic fleld at the nor- 
mal/superconducting interface, 

(V X h) X n|i = -Dnvnh,, (1.2) 

with Vn the interface velocity normal to the interface (with n ■ hj = 0); and a modifled 
Gibbs-Thomson boundary condition for the magnetic fleld at the interface. 
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where H^ is the thermodynamic critical field for the superconductor, Ons is the surface 
tension of the normal/superconducting interface, /C is the curvature of the interface, and 
r^^ is a kinetic coefficient for motion of the interface (here we ignore thermal fiuctuations 
1^). It is convenient to introduce the dimensionless surface tension a^s and kinetic coefficient 
r~^, which are defined through 

where A is the magnetic penetration depth, k is the Ginzburg-Landau parameter (the ratio 
of the penetration depth to the coherence length ^), m is the mass of a Cooper pair, and 
7 is a dimensionless order parameter relaxation time. The dimensionless surface tension 
and the kinetic coefficient can be expressed in terms of the solutions to the one dimensional 
equilibrium Ginzburg-Landau equations, which are 

1 



K,"^ 



F" -Q^F + F-F^ = 0, (1.5) 



Q" - F'Q = 0, (1.6) 

where F{x) is the dimensionless order parameter amplitude and Q{x) is the dimensionless 
magnetic vector potential [the magnetic field is H{x) = Q'{x)]. For an interface between 
the normal and superconducting phases, we have in the superconducting phase x —* — cxd, 
F(x) -^ 1 and Q{x) -^ 0, while in the normal phase x -^ oo, F(x) -^ and Q{x) ~ xj\f2. 
The surface tension and kinetic coefficient are then [H,^: 

ans = ^[/l(«:)-/2(/t)], (1.7) 

f-^ = h{K) - '^^^^^^h{K\ (1.8) 
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where the integrals I\ and I2 are defined by 
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h{K) = 2 dx {F'Y = 2k' / dx F'-F^- F'Q' , (1.9) 



/oo 1 /■oo 

dx -^Q'-{Q'f =2k' dxF'^Q', (1.10) 

-co y2 J—oo 

where second set of expressions are obtained by integrating by parts and using Eqs. ( |1.5|) 
and ( |1.6|) . The Ginzburg-Landau equations have analytical solutions only in certain limiting 
cases, to be discussed below. Somewhat surprisingly, there appear to be few numerical 
calculations of the surface tension, even though the fundamental importance of this quantity 
in distinguishing type-I and type-II superconductors was recognized by Ginzburg and Landau 
in 1950 |,|^. 

In order to complete the derivation of the free boundary model from the TDGL equations, 
in this paper we solve the equilibrium Ginzburg-Landau equations numerically for a wide 
range of k values, and use the solutions to calculate the surface tension and kinetic coefficient. 
Our paper is organized as follows. In Sec. II we review our numerical methods for solving 
the Ginzburg-Landau equations. In Sec. Ill we discuss the surface tension, and derive the 
asymptotic form of the surface tension for small-K, which agrees well with the numerical 
results for a large range of k values. Our results also show that an asymptotic expansion 
for the surface tension near k = 1/v^, which was derived in Appendix B of Ref. P], has 
a larger range of validity than expected. In Sec. IV we discuss our results for the kinetic 
coefficient. The kinetic coefficient can be either positive or negative depending upon the 
ratio of relaxation times which appear in the TDGL equations 0,0; from our results we are 
able to determine the range of parameter values which result in a positive kinetic coefficient. 
Sec. V is a discussion section in which we briefly summarize our results. 

II. NUMERICAL METHODS 



Numerical computation is necessary to solve Eqs. (|1.5| ) and (|1.6| ) in general. The method 
chosen was standard relaxation using Newton's method |^. Since the relative size of the 



solution scales as 1/k we rescaled the GL equations by substituting x' = kx which results 
in the rescaled equations 



Q^F + F-F^ 



0, 



(2.1) 



k'Q" - F'Q = 0. 



(2.2) 



The rescaled differential equations can then be written as the following set of first order 
finite difference equations: 



Ay2,fc/A4 
A|/3,fc/Ax'fc = yi^fc 

A2/4,fc/A4 



— 1/3, fej 



(2.3) 



with yi = F, ys = Q, ys = F' , y^ = Q', Ayn,k = yn,k - yn,k-i, Ax'^ = x'l, - x'^.^, and 
Fn fc = I iyn,k + yn,k-i)- Each of the four equations is to be solved for /c = 2 . . . M, with M 
the number of mesh points. We also have four boundary conditions as follows: 

2/1,1 = 1, 
2/4,1 = 0, 

yi,M = 0, 

1/4,M = 1/(k\/2), (2.4) 

giving a total of 4M equations of the 4M yn,kS. If we move the right hand side of Eqs. 
( p.3| ) to the left hand side and multiply by Ax'^ we are then left with a set of homogeneous 
equations. Labeling these equations En^k for n = 1 . . .4, k = 2 . . . M the problem is now to 
solve En^k = and Eq. ( |2.4| ) simultaneously. Given an initial guess y^ we can improve the 
solution using the expansion 

/ . ^ — ^yn,k-l + / „ 

^^l Oyn,k-l n=l ^yn,k 



4 op 4 o-p 

Efc(yfc + Ayfc, y^-i + Ayfc_i) ^ Efe(yfc, yfc_i) + ^ ^^_ ' Ayn,k-i + 51 ^::r~Ay„,fc, (2.5) 



where we want the left hand side to equal zero. This gives a system of linear equations to 
solve for the Ay^'s. We then add the Ay^ to the y^ to obtain a closer approximation. This 
process is repeated until the maximum value of \En,k\ is less than 10~^. 

The mesh of points x'^ are chosen at the start of the relaxation. We want the range 
to be large enough so that at the endpoints the functions are already close to their values 
at infinity. By trying different values for the endpoints we have found that at x' = ±100 
the functions are all sufficiently close to their boundary values at infinity that imposing the 



conditions ( p. 41) here do not significantly affect the solutions. We also want the mesh to 
be fine enough to accurately pick up the detail of rapidly varying areas. Knowing that the 
solutions all show the most change in a small area (taken to be around 0) and are relatively 
linear outside, we have manufactured a grid with more closely spaced points in the center 
than the edges using a total of 2001 points. The spacing was chosen so as to make the 
change in yi (i.e., F) from one mesh point to the other roughly constant. 

Convergence for this algorithm depends greatly upon the initial guess. For 0.1 < k < 10 
the solution is obtained within about 30 iterations from our initial guess. Since we are 
interested in finding solutions for a large range of k, it is advantageous to use the previous 
solution to start a new solution changing k slightly each time. In going from k =0.1 to 0.001 
by 0.001 each subsequent solution was obtained in only three iterations. 

The results of our computations for large and small k, are shown in Figs. (HI) and @. 
In Fig. (IID K = 10^^, and we see that the magnetic field is essentially a step-function, 
as suggested by the analysis in the next section. In this case the field only penetrates a 
short distance into the superconducting region, and therefore the full positive energy of fiux 
expulsion is obtained, resulting in a positive surface tension (type-I superconductor). In Fig. 
(0) K = 10; we see that the magnetic field penetrates far into the superconducting region, 
so that the positive energy of fiux expulsion is reduced. However, the negative condensation 
energy is very large here, resulting in a net negative surface tension (type-II superconductor). 
From the numerical solutions we computed the surface tension and kinetic coefficient using 
Eqs. ( pTTf) , ( |1.8|) , (|1.9|) , and (|1.10| ). Since our mesh spacing is already adapted to rapidly 
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varying parts of the solution, we used a basic trapezoidal rule to calculate the integrals 
Ji and l2- These results are tabulated in the Table. As expected, the surface tension 
passes through zero at k = 1/v^, separating type-I from type-II superconductors. In the 
following sections we will compare these numerical results against asymptotic solutions of 
the Ginzburg-Landau equations. 

III. SURFACE TENSION 

The surface tension is the excess free energy per unit area due to the presence of the 
interface. As shown by Ginzburg and Landau [Q (see also Ref. P), for k ^ I/v^j ^ns = 
2-\/2/3fi; + 0(fi:^^/^). Unfortunately, we have found that the lowest order expansion provides 
a very poor approximation to the surface tension, except for very small values of k (less than 
10^^); this fact was also noted by Ginzburg and Landau 0. Therefore, in this section we will 
first generalize their result somewhat by calculating the next order term in the expansion. 

In the small-K limit it is convenient to work with the rescaled Ginzburg-Landau equations, 
Eqs. ( p.l| ) and (|2.2| ). The lowest order approximation is obtained by setting the first term 



in the second Ginzburg-Landau equation, Eq. ( |2.2| ), equal to zero so that F^Q = 0. In the 
superconducting phase Q = with F in the superconducting phase determined by Eq. ( p.l| ) 
with Q = 0] the solution to this equation is F{x') = — tanh(x'/v2), for x' > 0. When 
this solution is substituted into the expression for the surface tension it produces the lowest 
order expansion derived by Ginzburg and Landau |^. To calculate the next order term, we 
need to take this expression for the order parameter and substitute it back into Eq. ( p. 21) , 
and then solve for Q. For x' > (the normal phase), we have Q'^ = 0, which integrates to 
Q>{x') = x' I K\pl + Ci, with C\ a constant to be determined by matching onto the x' < 
solution. For x' < (the superconducting phase), the vector potential satisfies 

K^Q"^ - tanh^(xVV2)g< = 0. (3.1) 

Although this equation does not appear to have an explicit analytical solution, for /t <^ 1 
we can use the WKB method to obtain an asymptotic solution. Some care is necessary as 



this equation has a second order turning point at x' = 0. The uniformly vahd asymptotic 
solution (i.e., a solution valid both near and away from the turning point) is [O] 



Q<{X') = C2 



25/8r(3/4) 1 



TT 



K 



1/4 



In cosh (x'/ a/2) 



1/2 



i^-1/4 



V2 



K 



lncosh(x'/v2 



(3.2) 



tanh(x'/y2)_ 

where C2 is a second constant of integration and K_iii^{z) is the modified Bessel function of 
order —1/4. The integration constants are determined by matching the solutions and their 
derivatives at x' = 0, with the result 



TT 



1 





^i - 


23/4r(3/4)2 0^' 




so that 




Q<{x') = 


1 

2V8r(3/4)fi;3/4 


"lncosh(x7y2)' 
_-tanh(x7y2)_ 


1/2 


^^lncosh(x7^) 



(3.3) 



(3.4) 

To calculate the surface tension, we substitute our solution into our expression for the 
surface tension, Eq. (|1.7| ). For I\ we have 



IAk) 





2k I dx' 

-00 



2^2 



F'Q 



2V2k 23/4 



TT 



and for lo 



/2(K)=2fi:/" dx'F^Q 



23/4 



TT 



81(3/4)2 



,3/2 



.3/2 



(3.5) 



81(3/4)2 

Therefore, from Eq. ( |1.7D we find that the surface tension in the small-K limit is 

2V2I 2^/^71 1 



(3.6) 



(Jn 



0(1) 



(3.7) 



3 K 4r(3/4)2y^ 

The first term in the expansion was previously obtained by Ginzburg and Landau 0, and 
the second term is the new result. This calculation can also be formulated as a variational 
calculation, with F{x') = — tanh(x'/^t,v^) a trial solution for the order parameter; the 
WKB calculation may be repeated, and the resulting solution used to calculate the surface 
tension as a function of the variational parameter C,v Minimizing this expression and taking 
the small-K limit, we obtain Eq. (|3.7| ) [jri|. 
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For K ^ '\-l \f^ the second derivative term in the first Ginzburg-Landau equation, Eq. 
( |1.5D , may be neglected, and the resulting algebraic equation solved for F as a function of Q. 
This expression is then substituted into the second Ginzburg-Landau equation, Eq. ( p.. 61 ), 
and the resulting nonlinear differential equation may also be solved ||^. The surface tension 
in this limit is o"ns = — 4(-\/2 — l)/3. 

The surface tension is zero at k = l/-\/2 0; at this point the Ginzburg-Landau equations 
become integrable [§]. The solutions may be used to carry out a local analysis of the surface 
tension about k = 1/v^ [Q, with the result a^s = 0.388(1/2^^ — 1). 

Summarizing, we have 

0.943/t-i - 0.880/t-i/^ K < 1; 

0.388(l/2fi:2 -1), ^ ^ 1/v^; (3.8) 

-0.552, /t > 1. 

In Fig. @ the numerical results for the surface tension are compared to the asymptotic 
expressions for k <^ 1. The asymptotic result is accurate for k < 0.2; the k~^/^ correction in 
Eq. (|3.7|) is important for values of k which are greater that 10~^. In Fig. (^) we compare 
the numerical results against the asymptotic expansion derived in Ref. for k ^ l/V^. 
The asymptotic expansion is reasonably accurate for 0.5 < n < LO. If we identify the small 
parameter in this expansion to be e = l/(2/i;^) — 1, then this would imply that the expansion 
is accurate for values of e as large as e = 1, a somewhat surprising result. In this figure we 
also see that the surface tension changes sign at k = l/-\/2, as expected. Fig. (^) shows the 
surface tension in the range 1.0 < k < 10.0; for large k we see that the surface tension is 
approaching the limiting value of —0.55. 

IV. KINETIC COEFFICIENT 

The kinetic coefficient T^^ is a function of k, which is a ratio of length scales, as well 
as 27rha^"^ / rwy , which is the ratio of the diffusion constant for the order parameter, D^ = 
h/2m'j, to the diffusion constant for the magnetic field, Dh = l/4:7ia^'^\ If this latter ratio 

9 



is sufficiently large, the kinetic coefficient may also change sign (resulting in some sort of 
dynamic instability). By setting T~^ = in Eq. (|1.8|) , we obtain an expression for the 
neutral stability curve : 



m'y 



27r^o-(") 



(4.1) 



neutral -'ll^'^jj 

The numerical result is plotted in Fig. (|^). 

In the limit of small-K we may use the previously derived expansions for /i and I2 in 
Eqs. (|3.5| ) and (|3.6| ) to obtain 



By setting T~^ = 0, we see that for small-K the stability curve should behave as k^^'^ for 
small-K, which is conffimed by the numerical results shown in Fig. (0). 



V. DISCUSSION 

We have investigated in some detail the behavior of the solutions to the one-dimensional 
Ginzburg-Landau equations using both numerical methods and asymptotic expansions. The 
numerical results for the surface tension of the normal-superconducting interface agree well 
with the small-K asymptotic expansion developed in this paper. In addition, a recently 
developed asymptotic expansion of the surface tension for values of k near l/-\/2 |^ agrees 
with the numerical results over a surprisingly large range of k values. We have also calculated 
the neutral stability curve for the kinetic coefficient, which will be important in studies of 
the dynamics of the normal-superconducting interface 
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FIGURES 
FIG. 1. Magnitude of the order parameter F and the magnetic field H for k = 10^'^, which 
corresponds to a type-I superconductor. Lengths are in units of the penetration depth. 

FIG. 2. Magnitude of the order parameter F and the magnetic field H for k = 10, which 
corresponds to a type-II superconductor. Lengths are in units of the penetration depth. 

FIG. 3. Dimensionless surface tension a^s as a function of the Ginzburg-Landau parameter k 
for 10^^ < K < 0.3. The solid line is the numerical result, and the dashed line is the asymptotic 



expansion for k <C 1 given in Eq. (3. 



FIG. 4. Dimensionless surface tension cjns as a function of the Ginzburg-Landau parameter k 
for 0.3 < K < 1.0. The solid line is the numerical result, and the dashed line is the asymptotic 
expansion about k = l/\/2 given in Eq. (|3.8|). 

FIG. 5. Dimensionless surface tension (Jns as a function of the Ginzburg-Landau parameter k 
for 1.0 < K < 10.0. The solid line is the numerical result, and the dashed line is the limiting value 
for K ^ 1 given in Eq. ( |3.8| ). 



FIG. 6. Stability diagram for the kinetic coefficient, determined by Eq. ( [4.1| ) in the text. The 
y-axis is m^/2iTha^'^\ the inverse of the dimensionless conductivity. For parameters in the region 
above the line the kinetic coefficient is positive, while the region below the line corresponds to a 
negative kinetic coefficient. 
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TABLES 



TABLE L Representative numerical results from the solution of the Ginzburg-Landau equations. 



K 


/i(k) 


^2W 


0"ns 


■p— 1 a 


0.001 


0.000926 


0.0000161 


910 


0.000910 


0.01 


0.00891 


0.000516 


84.0 


0.00840 


0.05 


0.0414 


0.00586 


14.2 


0.0356 


0.1 


0.0782 


0.0169 


6.13 


0.0613 


0.2 


0.144 


0.0495 


2.36 


0.0942 


0.3 


0.202 


0.0943 


1.19 


0.107 


0.4 


0.254 


0.150 


0.648 


0.104 


0.5 


0.301 


0.217 


0.338 


0.0845 


0.6 


0.345 


0.294 


0.142 


0.0511 


1/^2 


0.388 


0.388 


0.000219 


0.000109 


1.0 


0.490 


0.706 


-0.216 


-0.216 


5.5 


1.16 


17.2 


-0.529 


-16.0 


10.0 


1.43 


55.8 


-0.544 


-54.4 



We have chosen {2TTfia^"'' / rwy) = 1 for the purposes of illustration. 
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